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Abstract 

We discuss the relevance of higher order cumulants of net baryon number 
fluctuations for the analysis of freeze-out and critical conditions in heavy ion col- 
lisions at LHC and RHIC. Using properties of 0(4) scaling functions, we discuss 
the generic structure of these higher cumulants at vanishing baryon chemical po- 
tential and apply chiral model calculations to explore their properties at non-zero 
baryon chemical potential. We show that the ratios of the sixth to second and 
eighth to second order cumulants of the net baryon number fluctuations change 
rapidly in the transition region of the QCD phase diagram. Already at vanishing 
baryon chemical potential they deviate considerably from the predictions of the 
hadron resonance gas model which reproduce the second to fourth order cumulants 
of the net proton number fluctuations at RHIC. We point out that the sixth order 
cumulants of baryon number and electric charge fluctuations remain negative at 
the chiral transition temperature. Thus, they offer the possibility to probe the 
proximity of the chemical freeze-out to the crossover line. 
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1 Introduction 



Strongly interacting matter at high temperature or large net baryon number 
density is expected to undergo a rapid transition from a phase with hadrons 
as dominant degrees of freedom to a phase where partonic degrees of freedom 
prevail. At vanishing baryon chemical potential (fis = 0), this transition is a 
true second order phase transition only in the limit of vanishing light quark 
masses. For \ir > 0, however, a second order phase transition point, the 
so-called chiral critical point, may exist also for physical values of the quark 
masses. A large experimental as well as theoretical effort is put into the ex- 
ploration of the QCD phase transitions and the development of appropriate 
tools and observables that can provide a univocal signal for the existence of 
the phase transitions and their universal properties. The analysis of fluc- 
tuations of various physical observables p], in particular of the net baryon 
number, may serve this purpose [2j [3] . Theoretically, the properties of the 
corresponding susceptibilities are well understood at high temperature (T) 
and small values of the baryon chemical potential (/i_b). In this regime, they 
are suitable observables for localizing the phase boundary in the hb-T plane. 

Critical behavior is signaled by long range correlations and increased fluc- 
tuations, owing to the appearance of massless modes at a second order phase 
transition. Fluctuations of baryon number and electric charge have been 
shown to be sensitive indicators for such critical behavior [I]. In the ex- 
ploration of the QCD phase diagram at non-zero temperature and baryon 
chemical potential, higher order cumulants of baryon number fluctuations 
play a particularly important role. They diverge on the chiral phase transi- 
tion line T c {^Bi m q = 0) as well as at the elusive chiral critical point [5]. 

In heavy ion experiments, a lot of information has been collected on 
particle yields in a wide range of beam energies [6]. The particle multiplici- 
ties are well described in a thermal model using the partition function of a 
hadron resonance gas (HRG) [7j . Ratios of particle yields at a given beam 
energy can be characterized by a few thermal parameters, e.g. temperature 
and chemical potentials for baryon number, electric charge and strangeness. 
These parameters define the freeze-out conditions, i.e. the thermal param- 
eters corresponding to the last interaction of the hadrons participating in 
the collective expansion and cooling of the hot and dense matter formed in 
a heavy ion collision. Data obtained at small values of the baryon chemical 
potential suggest that the freeze-out curve Tf([i B ) is close to the expected 
QCD phase boundary. In particular, at /!# = the chemical freeze-out seems 



2 



to occur at or very near the QCD transition region for a physical quark mass 
spectrum [8]. At larger values of /j-b/T, however, there is a discrepancy be- 
tween the slope of the freeze-out curve and current lattice QCD results on 
the curvature of the chiral phase transition line [9] . 

The HRG model, which is based on the observed hadron spectrum, does 
not exhibit critical behavior nor does it reflect the sudden change of degrees 
of freedom in the transition to the partonic phase of QCD. In the chiral limit, 
close to the phase transition line T c (/i£,m g = 0) fluctuations of e.g. the net 
baryon number density are expected to reflect the universal properties [10] 
of the 3-dimensional, 0(4) symmetric spin model [11]. The 0(4) scaling re- 
lations for cumulants of net baryon number fluctuations differ significantly 
from predictions based on the HRG model. Lattice calculations of cumulants 
of baryon number and electric charge fluctuations, performed in the transi- 
tion region at vanishing baryon chemical potential and non-zero quark mass, 
do indeed show that these cumulants differ qualitatively from those of the 
HRG model and reflect the basic feature expected close to the chiral phase 
transition [4]. This suggests that at physical values of the light and strange 
quark masses these cumulants are sensitive to the critical dynamics in the 
chiral limit and may be employed to characterize also the crossover transition 
in strongly interacting matter. 

Thus, also at vanishing baryon chemical potential, i.e. under the condi- 
tions approximately realized in the high energy runs at RHIC or LHC, the 
question arises to what extent a refined analysis of the freeze-out conditions 
can establish the existence of the chiral phase transition. At /xb/T ~ 0, 
net quark number fluctuations and their higher cumulants can be computed 
within the framework of lattice QCD [H [12]. Such calculations will eventu- 
ally provide a complete theoretical characterization of the thermal conditions 
in the crossover region. This may then be used to unravel the relation of the 
freeze-out conditions at RHIC and LHC energies to the pseudo-critical line in 
the QCD phase diagram, provided the system remains close to thermal equi- 
librium during freeze-out. At present, however, lattice calculations provide 
only limited information on cumulants up to eighth order. In particular, 
controlled predictions on their properties in the continuum limit are still 
lacking; their characteristic features are obtained on a qualitative level, but 
quantitative results are not yet available. Viable alternatives for discussing 
qualitative features of the net baryon number fluctuations is offered by 0(4) 
scaling theory and by chiral models, like e.g. the Nambu-Jona-Lasinio (NJL) 
model. In particular, the effective models have the advantage that they can 
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be extended to /xg > with minimal effort. On the other hand, a clear 
disadvantage of NJL-type chiral models is that they do not account for the 
potentially large contribution from resonances in the hadronic phase. 

In this paper we will discuss the robust features of cumulants of net baryon 
number fluctuations that can be extracted from considerations based on 0(4) 
universality, on existing lattice calculations and on model calculations. In the 
next section we discuss higher order cumulants at vanishing baryon chemical 
potential making use mainly of 0(4) universality. In Section 3 we extend 
these considerations to hb/T > using results from model calculations. In 
Section 4 we summarize the relevance of our findings for experimental studies 
of baryon number fluctuations at LHC and RHIC and in Section 5 we give 
our conclusions. 



2 Charge fluctuations at /x# = 

2.1 0(4) scaling functions and net baryon number fluc- 
tuations 

Close to the chiral limit and at temperatures near the chiral phase transition 
temperature T c , higher order derivatives of the free energy density (/) with 
respect to temperature or chemical potential are increasingly sensitive to the 
non-analytic (singular) part (/ s ). We may represent the free energy density 
in terms of the singular and regular contributions 

f(T, fig, m q ) = f s (T, fi q , m q ) + f r (T, (i q , m q ) . (1) 

In addition to the dependence on temperature T, we also introduce here an 
explicit dependence on the light quark chemical potential, \i q = /is/3, and 
the (degenerate) light quark masses m q = m u = m^. For simplicity we do not 
take the chemical potentials for electric charge and strangeness into account, 
nor do we introduce an explicit dependence on the strange quark mass. The 
singular part of the free energy may be written as 

f'P'J^'V = Ah 1+1 ' 5 f f {z) , z = t/h 1 ^ 5 , (2) 
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where j3 and S are critical exponents of the 3- dimensional 0(4) spin model [TTj 
and 

1 (T-T c (ii q \? 



Here T c is the critical temperature in the chiral limit and to, ho are non- 
universal scale parameters (as is T c ). We use the chiral transition temperature 
T c also to set the scale for the explicit symmetry breaking, introduced by 
the non-vanishing light quark massed^- The amplitude A is fixed by the 
relation of the scaling function for the free energy density, ff(z), to the 
more commonly used scaling function fc{z), which characterizes the scaling 
properties of the chiral condensate, or in general the order parameter (M) in 
0(4) symmetric models, M = h 1 ^ 8 fc(z), where 

fc(z) = ~( 1 + l) //W + j/M) ■ (4) 

The scaling function ff(z) and its derivatives f^ n \z) have recently been 
determined for n < 3 using high precision Monte Carlo simulations of the 
3-dimensional 0(4) spin model and the known asymptotic series expansions 
[T6] . We will use these results as a starting point for a discussion of the generic 
structure of higher order cumulants of the net baryon number fluctuations. 

Note that the reduced temperature t, introduced in Eq. (J3J), depends ex- 
plicitly on the quark chemical potential. The constant K q ~ 0.06, which 
controls the curvature of the chiral phase boundary for small values of /x g /T, 
was recently determined in a scaling analysis of (2+l)-flavor QCD |9]. A 
comparison with other lattice results for n q [T71 [18], suggests that this pa- 
rameter is only weakly dependent on the quark mass and the number of 
flavors. 

In this paper we focus on the properties of the net baryon number fluc- 
tuations. The corresponding cumulants are obtained from Eq. ([TJ by taking 
derivatives with respect to fi q = fi g /T, 

b 1 d n f/T* 



1 In Ref. [15] the strange quark mass was used to set the scale for the symmetry breaking 
term. 



5 



From Eq. (J2J) it is apparent that, in the vicinity of the critical temperature, 
the susceptibilities Xn show a strong dependence on the explicit symmetry 
breaking term, the quark mass, 



-(2K q ) n / 2 h( 2 - a - n / 2 ^ s ff /2] \z) , for fi q /T = 0, and n even 
-(2«,) n fh^-^l^ffiz) , for fi q /T > 



where we used the scaling relation 2 — a = (35(1 + 1/5). Because a is 
negative in the 3-dimensional 0(4) universality class (a = —0.2131(34) [TT]). 
the specific heat, or equivalently the fourth order cumulants of the net baryon 
number fluctuations, does not diverge at the chiral transition temperature, 
i.e. at z — 0, in the chiral limit. At [i q fT = the first divergent cumulant 
is obtained for n — 6, while at fi q /T > this happens already for n = 3. 
We note that the singular structure appearing in n-th order cumulants for 
\x q jT > is identical to that of (2n)-th order cumulants at [iqjT = 0, since 
the chemical potential enters quadratically in the reduced temperature t. 
This characteristic has been used in Refs. [131 E] to exploit properties of 
third order cumulants at non-zero baryon chemical potential as signatures 
for critical behavior. In this case, however, the singular contributions are 
suppressed by a factor (ji q /T) n / 2 relative to the sixth order cumulants at 
fi q /T = 0. Consequently, in the mean-field analysis of Ref. [13], it was found 
that qualitative changes of the third-order cumulant, e.g. a change of sign in 
xf , is found only at rather large values of the chemical potential, Hb/T > 4. 
Using Eq. (Q, one finds the leading singularity in the chiral limit, 

B f-(2K 9 ) n / 2 |t| 2 - Q - n / 2 /± n/2) , for fi q /T = 0, and n even 



(2«,)» )" |t| 2 " Q -"/i nJ , for fi q /T > 



where 



= hm \z\-^-^ff\z) . (8) 

The singular part of xf, which is proportional to fj 2 \ has the same sin- 
gular structure as the specific heat; it is proportional to the second derivative 
of the free energy with respect to temperature. Thus, using the convention 
of Ref. [19J, we may write Xa (t) m the chiral limit at ji q /T = 0, 

Xf(*)-Xr + — \t\~ a , (9) 

a 



6 



where A + is the amplitude above and A~ below the critical temperature. The 
amplitudes A ± are positive and the ratio A + /A" ~ 1.8. This implies that 
the cumulants Xn(t) are positive for all n > 4 and t < 0, while for t > they 
alternate in sign. At non-zero values of the quark mass, h > 0, we thus expect 
Xq to change sign in the transition region and xi to do so twice. For a given 
h > O^this is reflected in the ^-dependence of the scaling functions fi (z), as 
shownzl in Figs. [TJandEJ m fact, the temperature and quark mass dependence 
of the singular parts of the net baryon number fluctuations is directly related 
to the scaling functions fi n \z) of the 3-dimensional 0(4) model [H] . Thus, 
the generic structure of the fourth and sixth order cumulants can be obtained 
from the known 0(4) scaling functions, in the chiral limit as well as for non- 
zero values of the quark mass. 

In the chiral limit, the non- analytic contribution to xf vanishes at the 
chiral transition temperature, t = 0. Consequently, in the transition region, 
the regular terms dominate in xf ■ Nonetheless, the non-analytic term in xi 
varies rapidly with temperature, leading to a pronounced maximum in the 
transition region, observed in lattice as well as model calculations. 

The temperature at which xi changes sign is non-universal since it de- 
pends on the magnitude of the regular terms. However, in the scaling regime 
the location of the extrema and the corresponding amplitudes follow univer- 
sal scaling laws. Moreover, we note that the positions of the two extrema of 
fj 3 \z), z~ < and z + > 0, provide bounds on the critical temperature in 
the chiral limit. 

It is evident from Figs. [TJ and [2] that the 0(4) scaling functions show much 
more structure and a stronger quark mass dependence in the symmetric 
phase, t > 0, than in the broken phase, t < 0. In the latter case, the 
divergence at t = builds up much more slowly than on the high temperature 
side. The scaling function ff(z) changes sign close to z — 0, while fj (z) 
does so already at z~ < 0. Furthermore, we note that the position of the 
maximum of ff\z) is at z + ~ 1.45 [TS], which is close to the peak position 

2 From the result for fj {z) |16| we also constructed an estimate for the next derivative, 

/^(z), which required some smoothening of the interpolations that entered the determi- 

nation of (z). The resulting scaling function and the resulting quark mass dependence 

of xi (*) i s shown in Fig. [2j We want to use it here to point out the qualitative structure 
that does arise within the 0(4) universality class, but do not consider this figure as being 
correct on the quantitative level, i.e. as far as the accurate location of the minima and 
the height of peaks is concerned. 
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Figure 1: Scaling of the non-analytic contributions to xf (left) and Xq 
(right) arising from second and third derivatives of the singular part of the 
free energy. Shown are results for different values of the symmetry breaking 
parameter h$h = m q /T c ; ho and zq = h^ 136 /to are non-universal scale pa- 
rameters. Note that for hoh = 1 the abscissa is the scaling variable z. The 
corresponding curve thus directly shows the 0(4) scaling function. 
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Figure 2: Same as Fig. [TJbut for the non-analytic contributions to xi ■ 



S 



of the chiral susceptibility, z p ~ 1.33(5) [UJ. The location of these extrema 
define pseudo-critical temperatures T 6± , which converge to that of the second 
order chiral phase transition temperature T c in the chiral limit, 

r6±(m « ) - 1 + z - faY" , do) 



where 1//35 ~ 0.55. The proportionality constant zq is uniquely determined 
in the scaling regime of QCD with two light quarks [15] and is the same 
constant, which also controls the scaling of the pseudo-critical temperature 
T x determined from the position of the peak in the chiral susceptibility, 

With decreasing quark mass the minimum of xi decreases as, 

, -(l+a)//8<J / \ -0.66 



"X,&,min ^ ( rp J ( rp J • (1^) 

We also note that corrections to ximin a ^ non-zero fi q /T start at 0((fi q /T) 4 ). 
We thus expect that the basic structure of higher cumulants persists also for 
I^b/T > 0. 



2.2 Sixth order cumulant and the QCD transition tem- 
perature 

Based on the generic structure of the 0(4) scaling functions we thus can 
understand the basic features of the temperature dependence of higher cu- 
mulants of the baryon number fluctuations at vanishing baryon chemical 
potential. We focus on the properties of the sixth order cumulant, x B {T), 
or correspondingly the ratio of cumulants Rq 2 (T) = xi (T) /xi (T) . In the 
hadronic phase, Xq{T) first grows with increasing temperature. It exhibits a 
maximum in the hadronic phase, close to the transition region and then drops 
rapidly. In the entire high temperature regime, xi{T) and consequently 
Rq 2 (T) remain negative. Lattice calculations of xfCO [H [20] suggest that 
in QCD with physical quark masses, these basic features, which are due to 
the singular terms in Xni persist. In particular, Xg(T) < in the vicinity of 
the pseudo-critical temperature for chiral symmetry restoration. 
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Figure 3: Schematic plot of the temperature dependence of the ratio of the 
sixth and second order cumulants of the net baryon number fluctuations in 
units of the phase transition temperature T c in the chiral limit. The vertical 
lines show the chiral phase transition temperature and the pseudo-critical 
temperature T pc , corresponding to a peak in the chiral susceptibility in QCD 
with physical light quark masses. 



We thus expect the sixth order cumulant to be a very sensitive probe for 
the temperature at which the freeze-out of hadrons in heavy ion collisions 
occurs. In fact, it is conceivable that hadronic freeze-out occurs in a tem- 
perature regime above the QCD phase transition temperature in the chiral 
limit, i.e. it may occur close to the pseudo-critical temperature that signals 
the onset of chiral symmetry restoration. In that case we would expect to 
observe 

Rg 2 (T) = ^flL-i < at freeze-out at LHC and RHIC high energy runs. 

This would be in striking contrast to conventional HRG model calculations, 
which predict Rq 2 = 1. A similar conclusion can be drawn for the eighth 
order cumulant xi or > equivalently, the ratio R$ 2 . 

In more general terms we suggest that a determination of Rg 2 at large 
collision energies at RHIC or LHC will provide a characteristic signature 
for the location of the freeze-out temperature relative to the QCD phase 
transition. This is illustrated in Fig. [31 where we show a schematic plot 
of Rf 2 along with the critical temperature in the chiral limit (T c ) and the 
crossover temperature for chiral symmetry restoration for physical light quark 
masses (T pc ). 
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In the following we extend the above considerations to the case of non- 
vanishing baryon chemical potential. We do this in the framework of a chiral 
model, analyzing the critical behavior of the Polyakov loop extended quark 
meson (PQM) model in the functional renormalization group approach [14] . 
This allows us to compute the dependence of higher cumulants on fi q /T and 
T. We determine the line of minima of R§ 2 m the \i q -T plane and compare 
this with the chiral transition line as well as with the pseudo-critical line 
obtained at a physical pion mass. 

3 Higher cumulants of charge fluctuations at 
Vb/T > 

The universal features of the 3-dimensional 0(4) scaling functions discussed 
in the previous section also are reflected in effective chiral model calcula- 
tions where the Polyakov loop is coupled to the fermion sector, thus generat- 
ing many of the characteristics of the connnement-deconfinement transition 
[2H [221 [231 [21]. The Polyakov loop extended quark meson (PQM) model [25] 
is one variant. It shares with QCD a global 0(4) symmetry, which is spon- 
taneously broken at low and restored at high temperatures. Thus, in the 
light quark mass limit, this model reproduces the universal scaling functions, 
discussed in the previous section. We use the PQM model to implement 
the basic features of the transition from hadronic matter to a quark gluon 
plasma into the temperature dependence of the cumulants of net baryon 
number fluctuations and to analyze their properties also at [i q jT > 0. 
The Lagrangian of the PQM model [25] . 

L = q[i 1 ^-g(a + i 1 ,TK)]q+ l -{d,a) 2 + l -{d^) 2 

-U(a,7r) -U{i,f) , (13) 

involves interactions of mesons and gluons with fermionic fields. The coupling 
between the effective gluon field and quarks is implemented through the 
covariant derivative 

= - (14) 

where A^ = g s A a ^\ a /2 and the spatial components of the gluon field are 
neglected, i.e. A^ = S^qAq. 
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The purely mesonic potential of the model, 



U(a, 7f) = j {a 2 + t? 2 -v 2 ) 2 - ca, (15) 

contains a term linear in the a field, which explicitly breaks the chiral sym- 
metry and is used to reproduce the physical pion mass. 

The effective potential for the gluon field is expressed in terms of the 
thermal expectation values of the color trace of the Polyakov loop and its 
conjugate, 

U(£,t) ^(T)^_63 (£3 + r3) + 64 (r£)2 (lg) 



where 



and 



T 4 2 6 v ' 4 



^<Tr c L(x)>, r = ^(Tr c Lt(x)>, (17) 



L(x) =V exp 



drAJx, t) 



o 



(18) 



Here V stands for the path ordering, /3 = 1/T and A4 = i Aq. The potential 
U(£,£*) preserves the Z(3) symmetry of the gluonic sector of QCD. 

In the Appendix, we give further details on the choice of parameters 
for the Polyakov loop potential and present the relevant calculational steps 
within the functional renormalization group (FRG) approach [26 | [27 ] |28 | 129]. 
Within this formalism we compute the free energy density of the PQM model 

[12 ED], 

f PQM (£, t; T, aO = n(£, t; T, /i) = ~ In Z(£, £*; T, , (19) 

as a function of temperature, chemical potential as well as the Polyakov 
loop variables, £ and £*. The latter two are then fixed by the stationarity 
condition: 

^n{£X;T,iM) = o, An(f ) r;T )/1 ) = o. (20) 

The cumulants of the net baryon number fluctuations are obtained by taking 
suitable derivatives of the free energy, as specified in Eq. (J5]). In practice 
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Figure 4: The sixth and eighth order cumulants of the net baryon number 
fluctuations at fi q /T = in the PQM model. The temperature is given in 
units of the pseudo-critical temperature T pc (m n ) corresponding to a maxi- 
mum of the the chiral susceptibility. The shaded area indicates the chiral 
crossover region. 

these derivatives have been implemented directly into the analysis of the 
flow equations (see Appendix). 

In Fig. H]we show the sixth and eighth order cumulants of the net baryon 
number fluctuations computed at \i q jT = within the PQM model for phys- 
ical values of the pion mass. The basic features dictated by 0(4) symmetry 
restoration, as discussed in the previous sections, are readily identified in the 
figure. Moreover, the positions of the two extrema of xi correspond approx- 
imately to the zeros of xi ■ This confirms that in the transition region, two 
derivatives with respect to [i q jT are indeed equivalent to one derivative with 
respect to T. 

From these calculations, as well as from calculations of the lower order 
cumulants X2 an d Xa j we obtain the ratios R^m °f the n-th and m-th cu- 
mulants. Results obtained for fi q /T = and fi q /T > are shown in Figs. \5\ 
and El respectively. We note that these ratios approach unity at low tem- 
peratures, as it is the case also in the hadron resonance gas model. In the 
transition region, they reflect the expected 0(4) scaling properties; they have 
a shallow maximum close to the transition region before they drop sharply. 
In particular, they show pronounced minima with R^ 2 < in the vicinity 
of the chiral crossover temperature. The exact location of these minima and 
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Figure 5: The temperature dependence of the fourth, sixth and eighth order 
cumulants of the net baryon number fluctuations Xn relative to the second 
order one. The temperature is given in units of the chiral crossover temper- 
ature. The shaded area indicates the region of the chiral crossover transition 
at fJjq/T = 0. The calculations were done in the PQM model within the FRG 
approach. 



their depth is to some extent model dependent. However, we note, that in 
the transition region the second order cumulant used in these ratios for nor- 
malization is dominated by non-singular contributions which are positive. 
The minima in R^ 2 therefore mainly reflect the strong temperature depen- 
dence of higher cumulants Xn- We also note that these minima become more 
pronounced with increasing fi q /T. In fact, the structure of e.g. R$ 2 becomes 
similar to that of xi a ^ large fi q /T. This is easily understood in terms of the 
Taylor expansion of R$ 2 , where the dominant correction at non-zero fi q /T is 
due to xi , 

KMT) = < 2 (0) + i(^) 2 (^ 2 (0)-< 2 (0)^ 2 (0)) 

+0{{ N /Tf) . (21) 

This also makes it clear why for fi q /T > the location of the minimum 
of RQ 2 (fi q /T) is shifted to lower temperatures relative to that of the chi- 
ral crossover temperature. Similarly, at non-zero ft q /T, the ratio Rf 2 (fi q /T) 
shows more pronounced oscillations in the transition region, due to contri- 
butions from higher order cumulants, which oscillate more rapidly in the 
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Figure 6: Comparison of the temperature dependence of the ratios Xe /xi 
and xf/xf f° r various fi q /T corresponding to values at chemical freeze-out 
in heavy ion collisions at RHIC. The shaded area indicates the region of the 
chiral crossover transition at fi q /T = 0. 



transition region. The amplitude of the maximum at high temperatures be- 
comes compatible in magnitude with that of the minimum. 

4 Freeze-out and the QCD transition 

The analysis of universal scaling functions that control the thermodynamics 
in the vicinity of a phase transition in the universality class of 3-dimensional 
0(4) symmetric theories (section 2), of model (section 3) and of lattice cal- 
culations [201 GS| su gg es t that at vanishing baryon chemical potential the 
sixth order cumulant of the net baryon number fluctuations is negative in 
the entire high temperature phase. In fact, the 0(4) scaling functions for 
higher cumulants turn negative already in the vicinity of the chiral (m q = 0) 
critical temperature, i.e. below the crossover temperature (m q > 0), which 
is relevant for the transition at non-zero values of the light quark masses. 
The regular terms in the QCD free energy may shift the onset of the neg- 
ative regime to higher temperatures. However, model [32] as well as lattice 
[20| [33] calculations suggest that the regime of negative sixth order cumulants 
starts below but close to the QCD crossover transition temperature. 

At non-zero baryon chemical potential, the temperature interval of nega- 



15 



H 
H 



1. 

£0. 


0, 




1 1 1 


1 1 1 1 


1 1 1 1 1 
ch. cross. _ 






minJCg/Xj 




~\ 






\ 








s. \ 

\ N 
\ \ 






\ \ 
\\ 






\\ " 



0.2 0.4 0.6 0. 

u IT 

rq pc 



1 1.2 



1.2 
1.1 



ch. cross. 

B, B 

mm x 8 /X 2 




Figure 7: The chiral crossover line [dashed line] and the first minima in \& 
(left) and xf (right) [solid line]. The bands show the parameter range for 
which x& an d xi > respectively, are negative in the neighborhood of these 
minima. 



tive xi > or equivalently R^i^q/T), shrinks and follows the crossover transi- 
tion line. This is illustrated in Fig. which shows the temperature interval, 
closest to the hadronic phase, where the sixth and eighth order cumulants 
of the net baryon number fluctuations are negative, as obtained in the FRG 
approach to the PQM model. It is evident that the sixth order cumulant xi 
is negative in a wide range of temperatures which extends into the symme- 
try broken phase. This is even more the case for the eighth order cumulant 
as expected from the structure of the corresponding 0(4) scaling function. 
Except for a small range of chemical potential values close to fi q /T = 0, the 
eighth order cumulant is, however, positive again on the crossover line. 



5 Discussion and Conclusions 

We have shown that higher order cumulants of the net baryon number fluc- 
tuations are sensitive probes for the analysis of freeze-out conditions in heavy 
ion collisions and may allow to clarify their relation to the QCD phase transi- 
tion. This is the case at LHC energies as well as at the entire regime of beam 
energies covered by the low energy run at RHIC. If in heavy ion collisions, 
particles are produced from a thermalized system, the analysis of higher cu- 
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mulants of the net baryon number fluctuations does provide constraints on 
the location of the freeze-out temperature relative to the chiral transition 
temperature. We find that the most robust statements, which become rigor- 
ous in the chiral limit, can be made for the sixth order cumulants of charge 
fluctuations for small values of the quark chemical potential: 

If freeze-out occurs close to the chiral crossover temperature the sixth 
order cumulant of the net baryon number fluctuations will be negative at LHC 
energies as well as for RHIC beam energies ^/snn ^ 60 Ge V, corresponding 
to /ib/T<0.5. This is in contrast to hadron resonance gas model calculations 
which yield a positive sixth order cumulant. 

We also note that the basic features discussed here for net baryon number 
fluctuations also carry over one-to-one to electric charge fluctuations. In the 
chiral limit, the most singular component in the cumulants of electric charge 
fluctuations is proportional to the singular part of the corresponding cumu- 
lant of net baryon number fluctuations. In fact, lattice [33] and model [M] 
calculations of the sixth order cumulant of electric charge fluctuations at 
H q /T = show that this cumulant is negative in the high temperature phase 
of QCD. The ratio of the sixth and second order cumulants, Xe /X2 rapidly 
drops in the transition region from the HRG value (Xe/X2 )hrg — 10 to 
about —4 at the chiral crossover temperature. Also this ratio is therefore a 
sensitive probe of the conditions at freeze-out and their relation to the critical 
behavior in strongly interacting matter. 

We finally comment on the fourth order cumulants, in particular on the 
ratio Xa 1X2 recently measured by STAR [35] as well as the corresponding 
ratio for cumulants for electric charge fluctuations, x^/xi- The former is 
consistent with the HRG value (xf /xi )hrg — 1 [251 EI]- ^ freeze-out occurs 
in the hadronic phase, one expects to find also large values for ratios of 
cumulants of electric charge fluctuations, {Xi/xi )hrg — 1.8 [36]. On the 
other hand, if freeze-out occurs at the crossover temperature, also these ratios 
will deviate from the HRG values. Lattice calculations suggests that both 
ratios drop rapidly in the transition region, but stay positive also in the high 
temperature phase [55] . 

In Table [1] we summarize two different freeze-out scenarios characterized 
by various ratios of the cumulants of the net baryon number and electric 
charge fluctuations, respectively. These scenarios assume that in heavy ion 
collisions the freeze-out of hadrons occurs from a thermalized system charac- 
terized by the phenomenologically determined freeze-out curve Ty(/x^). Ta- 
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freeze-out conditions 


x?/xf 


xf/xf 


x 4 y /x^ 


xt/xt 


HRG 


1 


1 


~ 2 


~ 10 


QCD: T freeze /T pc <0.9 


>1 


>1 


~ 2 


~ 10 


QCD: Tf reeze /T pc ~ 1 


~ 0.5 


< 


~ 1 


< 



Table 1: Values for ratios of cumulants of net baryon number (B) and elec- 
tric charge (Q) fluctuations for the case that freeze-out appears well in the 
hadronic phase (third row) or in the vicinity of the chiral crossover tem- 
perature (fourth row). We give results based on current lattice calculations 
[20], [33] and on the calculations presented here. In the second row we give re- 
sults of a HRG model calculation [33]. We also note that unlike the cumulants 
of net baryon number fluctuations the ratios of cumulants of electric charge 
fluctuations vary somewhat as a function of the baryon chemical potential 
along the freeze-out line. 

ble[T] shows, that ratios of cumulants of charge fluctuations are very sensitive 
to possible differences in freeze-out and crossover temperature. Finally, we 
note that these results do not account for possible finite volume effects nor 
for possible effects of the evolution from chemical towards thermal freeze-out 
in heavy ion collisions. 
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Appendix 

Here we summarize the basic steps involved in the calculation of the free 
energy of the PQM model within the functional renormalization group (FRG) 
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approach. We also describe the methods used to obtain the FRG results for 
the cumulants of the net baryon number fluctuations. Further details on 
these calculations are given in Refs. [301 E2 ED] . 

In the PQM model Lagrangian introduced in Eqs. (I13H18I) . we need to 
specify the parameters used in the Polyakov loop potential ffTB"]) . These pa- 
rameters were chosen to reproduce the equation of state of pure SU (3) lattice 
gauge theory, 



b 2 (T) 



a + a\ 



T 



+ a 2 



To 
T 



+ a 3 



To 
T 



2.625, a 3 = -7.44, b 3 = 0.75, b A 



(22) 
7.5 and 



with a = 6.75, ai = —1.95, a 2 
T = 270 MeV. 

The FRG flow equation for the thermodynamic potential in general con- 
tains mesons, quarks and Polyakov loops as dynamical fields. However, in 
the current calculation, we treat the Polyakov loop as a background field 
which is introduced self-consistently on the mean-field level. Following pre- 
vious work pH [30] , we formulate the flow equation for the scale-dependent 
grand canonical thermodynamic potential density for quarks and mesons 



d k Q k (£,t;T, N ) 



k 4 



12tt 2 EL 



2n B {E^T) 



1 



2n B {E a ,T) 



AN r N 



f 



Here 



E„ 



n B (E n y,T) 



N(£,t;T,ii q )-N(e,e*;T,v q ) 



.(23) 



1 



exp(E 7rjCT /T) - 1' 

is the bosonic distribution function, with the pion and sigma energies, 



EL 



k 2 + n k , E a = dk 2 + n k + 2 P n 



-a 



(24) 



(25) 



where Q k = &>k + ca and primes denote the derivatives with respect to 
p = a 2 /2. The momentum distributions of quarks N and antiquarks N 
are modified owing to the coupling to gluons, 

1 + 2ng£* + n 2 J 



N(l,t\T,n q ) 
N{£,t;T, N ) 



1 + 3n 2 gl 



N{£*,£;T, 



(26) 
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where E q = \Jk 2 + 2g 2 p is the quark energy and riq(E q ; T, p q ) = exp[(3(E q — 

We solve the flow equation for Qk{T,p q ) numerically by expanding it in 
powers of field variable, p, around the scale dependent potential minimum 



M 



n k (T,p q ) = 



^m,k(T, p q ) 



m=0 



ml 



(p - PkY 



(27) 



The location of the scale dependent minimum of thermodynamic potential 
flk{T,p q ) is determined by the stationarity condition 



da 



dVt k 



a=<T k 



da 



c = 0. 



a=a k 



Truncating the expansion in Eq. (|2"T|) at M — 3, Eq. 
equation for the Taylor coefficients and field variable, pk, as 



dkCtO,k 

d k Pk 

9kO-2,k 
9kCl3,k 



/ 2f>k 



dkpk + <9*A 



c/{2p k fl 2 + a 2 , k 
a 3,k dkpk + dk^l'k, 



d k £l' k , 



(28) 

yields the flow 

(29) 

(30) 

(31) 
(32) 



These differential equations are solved numerically with the initial cutoff 
A = 1.2 GeV (for details see Refs. [T^l 130]). The initial conditions for the 
flow are chosen to reproduce the pion mass in vacuum m n = 138 MeV, the 
pion decay constant /„- = 93 MeV, the sigma mass m a = 600 MeV and the 
constituent quark mass m q = 310 MeV at the scale k = 0. 

By construction, the solution of Eq. (1231) . Q k , computed at the mini- 
mum p = pk describes the scale-dependent thermodynamic potential den- 
sity, where the contribution of soft quark and meson modes with |g| < k is 
suppressed [3D]. Integrating the flow equation from k = A to k — >• 0, we ob- 
tain the thermodynamic potential which accounts for all momentum modes 
\q\ < A. The solution of equations ( 12911321) . yields the thermodynamic po- 
tential density for quarks and mesons, Q k ^o(£, £*'■, T, function of the 
Polyakov loop variables £ and £*. The full thermodynamic potential density 
fl(£, £*; T, p q ) in the PQM model, including quarks, mesons and gluons, is 
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obtained by adding the effective gluon potential to Q^oC^ T, fi g ), 

n(£, t; T, N ) = n k ^ (£, t; T, N ) + U(£, t), (33) 

where at a given temperature and chemical potential, the Polyakov loop 
variables are determined by the stationarity conditions: 

^n(*,r ; :r,//,) = o, An(^r ; r,/x,) = o. (34) 

The thermodynamic potential (!33l) does not contain the contributions of 
statistical modes with momenta larger than the cutoff A. In order to obtain 
the correct high-temperature behavior of the thermodynamics we supplement 
the FRG potential with the contribution of high- momentum states with \q\ > 
A. This contribution to the flow is approximated by that of massless quarks 
interacting with the Polyakov loops [TH [30] , 

d k Q k {T,fiq) = (35) 

l-N(£,t;T, N )-N{£,t;T,^ q ) . 

The complete thermodynamic potential of the PQM model is then ob- 
tained by integrating Eq. fl35|) from k = oo to k = A, where we switch to the 
PQM flow equation (1231 . The thermodynamic potential thus obtained can 
be used to explore the properties of different cumulants of the net baryon 
number fluctuations. 

The flow equation yields the thermodynamic potential as a function of 
temperature and chemical potential. Different cumulants of the net baryon 
number fluctuations x^, Eq- ©, can then in principle be obtained by explicit 
numerical differentiation. However, owing to the numerical errors accumu- 
lated in the solution of the flow equations, this method does not provide 
reliable results. Therefore, we compute the higher cumulants by analytically 
differentiating the flow equations ( |29H32|) and (|35|) with respect to fi q . This 
defines flow equations for \i an d X2 ■> which are solved using the solution of 
Eqs. ( 1291I321 as input. Subsequently, the higher cumulants of the net baryon 
number fluctuations up to eighth order are computed by numerical differen- 
tiations of X2 ■ 
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